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ABSTRACT 

This  paper  presents  new  capabilities  of  the  virtual-human 
Santos™  introduced  last  year.  Santos™  is  an  avatar  that 
has  extensive  modeling  and  simulation  features.  It  is  a 
digital  human  model  with  over  100  degrees-of-freedom 
(DOF),  where  the  hand  model  has  25  DOF,  direct 
optimization-based  method,  and  real-human  like 
appearance.  The  newly  developed  analysis  includes  (1) 
a  25-DOF  hand  model  that  is  the  first  step  to  study  hand 
grasping;  (2)  posture  prediction  advances  such  as 
multiple  end-effectors  (two  arms,  two  arms  +  head  + 
legs),  real-time  inverse  kinematics  for  posture  prediction 
for  any  points,  vision  functionality;  (3)  dynamic  motion 
prediction  with  external  loads;  and  (4)  musculosteletal 
modeling  that  includes  determining  muscle  forces,  and 
muscle  stress.  With  these  newly  developed  capabilities 
Santos™  can  be  used  to  test  the  joystick  design,  study 
grasping,  facilitate  vehicle  interior  design,  test  visibility  for 
product  design,  predict  correct  dynamic  motion  or 
posture  subject  to  external  loads,  and  investigate  muscle 
forces,  and  muscle  stress.  Finally,  additional  ongoing 
projects  are  summarized. 

Keywords:  Virtual  humans,  direct  optimization-based 
approach,  human  performance  measures. 

INTRODUCTION 

The  development  of  realistic  virtual  (computer-based) 
humans  is  quickly  gaining  momentum  in  the  research 
arena  as  well  as  in  the  industry  and  military  sectors. 
There  are  several  commercial  products  available  such  as 
Ramsis®,  Jack®,  Safework®,  Realman®,  etc.  None  of 
them  has  been  developed  using  the  direct  optimization- 
based  approach.  In  last  year's  report  we  presented  the 
new  generation  of  virtual  humans  (Santos™).  This  paper 
will  introduce  newly  developed  capabilities  of  Santos™. 


There  are  essentially  three  stages  to  developing  virtual 
humans:  (1)  basic  human  modeling,  (2)  input 
functionality,  and  (3)  intelligent  reaction  to  input  (memory, 
reasoning,  etc.).  This  paper  and  last  year's  report 
address  developments  in  the  first  stage  of  virtual 
humans.  The  ultimate  goal  will  be  towards  the  second 
and  third  stages. 

This  paper  presents  the  newly  developed  capabilities  for 
Santos™  that  entails  hand  modeling;  posture  prediction 
with  multiple  end-effectors;  inverse  kinematics  (IK);  a 
vision  objective  function,  and  self-collision  avoidance; 
motion  prediction  with  external  loads;  and 
musculoskeletal  modeling. 

HAND  MODEL 

A  virtual  human  hand  is  an  essential  part  of  a  virtual 
human's  form,  function,  and  communication,  capable  of 
complex,  expressive  articulation.  There  are  considerable 
works  on  hand  modeling.  ElKoura  and  Singh  (2003) 
categorized  hand  modeling  as  engaging  several  areas  of 
research  such  as  anatomy,  robotics,  animation  industry, 
graphics  and  vision,  and  music.  Wan  et  al  (2004),  and 
Savescu  et  al  (2004)  also  developed  different  hand 
models.  The  hand's  degrees  of  freedom  (DOF)  are  from 
21  to  27  for  different  research  groups  in  the  literature. 
Our  objective  is  to  develop  a  realistic  hand  model  based 
on  the  work  of  previous  researchers. 

A  human  hand  has  27  bones  (Fig.  1).  The  symbolic 
Roman  numbers  are  assigned  for  fingers:  the  thumb  is  I, 
the  index  finger  is  II,  the  middle  finger  is  III,  the  ring  finger 
is  IV,  and  the  little  finger  is  V.  At  the  base  of  each  finger 
is  one  metacarpal  bone,  which  connects  to  the  wrist. 
Carpometacarpal  (CMC)  joints  connect  the  metacarpal 
bones  to  the  wrist.  The  CMC  joints  for  II  and  III  are  static: 
they  cannot  actually  rotate.  In  contrast,  the  CMC  joints  for 
I,  IV,  and  V  can  rotate  with  Flexion/Extension  (F/E)  and 
Abduction/Adduction  (Ab/Ad).  The  metacarpophalangeal 
(MCP)  joints  connect  the  metacarpal  bones  to  the 


the  metacarpal  bones  to  the  phalanx  bones.  These  MCP 
joints  can  also  move  with  2  DOF  (F/E  and  Ab/Ad).  An 
interphalangeal  (IP)  joint  connects  the  two  phalanx  bones 
in  the  thumb  and  has  one  rotational  degree  of  freedom. 
Proximal  Interphalangeal  (PIP)  joints  connect  the  two 
inner  (not  the  tip)  phalanxes  in  II,  III,  IV,  and  V.  These 
joints  have  one  rotational  DOF  in  the  F/E  direction.  The 
last  joint  is  the  distal  interphalangeal  (DIP)  joint.  It 
connects  the  fingertips,  the  final  phalanxes  in  each  of  the 
fingers  (II,  III,  IV,  and  V),  and  it  moves  in  F/E  direction. 


From  the  above  described  anatomy,  the  hand  can  be 
modeled  as  25  DOF  shown  in  Fig.  2.  All  of  joints  in  the 
hand  are  revolute  joints. 

Hand  bones  have  different  lengths  for  different  fingers. 
All  link  lengths  depend  on  two  parameters  HL  and  HB , 
where  HL  is  the  length  of  the  hand,  and  HB  is  the  width 
of  the  hand  (Pena  Pitarch  etal 2005). 

ADVANCES  IN  POSTURE  PREDICTION 

There  are  several  approaches  for  posture  prediction  such 
as  empirical  approaches  (Beck  and  Chaffin  1992;  Berme 
etal,  1987;  Das  and  Behara  1998;  Faraway  et  al,  1999), 
inverse  kinematic  methods,  in  particular  pseudo-inverse 
methods  (J  ung  and  Choe  1996;  J  ung  et  al,  1992;  1995; 
Kee  et  al,  1994),  and  the  optimization-based  method 
(Yang  et  al,  2004).  It  is  important  to  consider  multiple 
end-effectors,  collision-free,  and  incorporated  vision 
objective  function  in  the  optimization-based  method  for 
posture  prediction. 

In  this  section,  we  first  provide  a  brief  overview  of  our 
approach  to  posture  prediction.  We  then  discuss  the 
basics  of  some  of  the  latest  developments  of  this  area 
with  respect  to  Santos™  Four  key  advancements  are 
presented:  1)  the  ability  to  predict  postures  with 
restrictions  on  any  number  of  end-effectors  on  the  body, 
2)  an  innovative  version  of  human  inverse  kinematics,  3) 
a  new  performance  measure  that  governs  predicted 
postures  based  on  one's  ability  to  see  an  object,  and  4) 
incorporation  of  self-avoidance,  a  capability  that 
essentially  makes  Santos™  aware  of  his  own  body. 


Fig.  1.  Right  hand  with  5  rays  (in  a  posterior  view) 


Fig.  2.  Dorsal  view  right  hand  model 
GENERAL  APPROACH 

Our  direct  optimization-based  approach  to  posture 
prediction  is  summarized  here,  while  Marler  et  al  ( 2005a, 
2005b)  provide  details.  This  approach  involves  finding 
the  joint  angles,  q,  in  order  to  minimize  one  or  more 
human  performance  measures  that  serve  as  the 
objective  functions  and  govern  posture.  The  problem  is 
constrained  by  limits  on  the  joint  angles  and  by  the 
requirement  that  the  avatar  must  contact  the  target  point. 
The  final  optimization  formulation  is  given  as  follows: 

Find:  qe  RD0F  (1) 

to  minimize:  Human  performance  measure(s) 

subject  to:  [x(q 

qf<qi<qf-i  =  \,2,...,DOF 


This  problem  is  solved  numerically  using  the  software 
SNOPT  (Gill  et  al,  2002),  which  uses  sequential 
quadratic  programming  (Arora,  2004).  Details 
concerning  the  different  components  of  the  optimization 
formulation  are  discussed  as  follows. 

The  design  variables  are  the  generalized  coordinates  qt , 

which  represent  the  components  of  the  vector  q.  These 
design  variables  have  units  of  radians.  The  first 
constraint  in  (1)  is  the  distance-constraint,  and  it  requires 
that  the  end-effector  contact  a  predetermined  target  point 
in  Cartesian  space,  where  e  is  a  small  positive  number 
that  approximates  zero.  In  addition,  each  joint  angle  is 
constrained  to  lie  between  lower  and  upper  limits, 
represented  by  qf  and  ^respectively.  These  limits 
may  be  changed  according  to  different  sets  of 
anthropologic  data. 

One  of  the  most  significant  advantages  of  this  approach 
is  that  it  is  completely  predictive;  it  does  not  depend  on 
predetermined  data  or  prerecorded  simulations.  It  thus 
affords  the  avatar  a  substantial  amount  of  autonomy.  In 
addition,  this  approach  is  relatively  fast,  providing  results 
in  real  time.  Finally,  various  capabilities  and  types  of 
functionality  can  be  incorporated  simply  by  modifying  or 
adding  different  performance  measures  and/or 
constraints.  For  instance,  multiple  end-effectors  can  be 
modeled  with  additional  constraints  (Farrell  et  al,  2005). 
Orientation  of  various  body  parts  can  be  specified  with 
additional  constraints  as  well  (Yang  et  al,  2005).  Ideas 
such  as  musculoskeletal  discomfort  can  be  modeled  as  a 
performance  measure  (Marler  et  al,  2005).  In  the 
following  sections,  we  discuss  how  multiple  end  effectors 
and  self-avoidance  can  also  be  incorporated  as 
constraints,  and  how  the  effects  of  vision  on  posture  can 
be  studied  using  a  new  performance  measure. 
Ultimately,  these  posture  prediction  capabilities  serve  as 
a  component  of  dynamic  motion  prediction,  which  is 
discussed  later  in  the  paper. 

MULTIPLE  END-EFFECTORS 

As  suggested  in  Yang  et  al  (2004),  the  DH-method  is 
used  to  determine  the  Cartesian  location  of  an  end- 
effector,  based  on  joint  angles.  This  method  can  be  used 
to  determine  the  position  of  any  point  on  the  body.  In 
fact,  a  user  can  specify  any  point(s)  in  space  relative  to 
any  single  local  coordinate  system(s).  Such  point(s)  can 
then  be  constrained  to  contact  any  target  point(s),  by 
adding  additional  distance  constraints  in  (1),  as 
demonstrated  in  Fig.  3.  This  capability  is  used,  for 
example,  when  Santos™  is  required  to  have  his  elbows 
remain  on  armrests  or  his  back  remain  against  the  back 
of  a  seat. 


Fig.  3.  Multiple  End-effectors 

ADVANCED  INVERSE  KINEMTICS  (IK) 

Virtual  humans  are  most  commonly  used  as  mannequins 
that  can  be  positioned  manually  in  a  particular  setting  or 
in  virtual  prototype.  It  is  crucial  that  Santos™  not  lose  this 
necessary  functionality.  To  be  sure,  the  predictive 
capabilities  of  Santos™  exceed  those  of  currently 
available  tools,  nonetheless,  a  user  must  be  able  to 
position  Santos™  easily.  We  have  developed  two 
versions  of  inverse  kinematics  (IK):  standard  IK  and 
advanced  IK.  Especially  in  advanced  IK,  we  incorporate 
joint  range  of  motion,  and  the  results  not  only  include 
joint  angles,  but  also  entail  the  values  of  human 
performance  measures  in  real-time. 

With  standard  IK  (used  in  Maya®),  the  user  is  able  to 
select  hot  points  on  Santos™  and  then  place  them 
wherever  is  necessary,  as  shown  in  Figs.  4  and  5.  The 
complete  avatar  then  moves  accordingly,  in  real  time. 
However,  standard  IK  does  not  incorporate  joint  ranges 
of  motion  and  sometimes  it  generates  odd  postures. 
Therefore,  we  extend  direct  optimization-based  posture 
prediction  to  advanced  IK. 

With  advanced  IK,  we  capitalize  on  posture  prediction 
capabilities  to  yield  an  exciting  new  tool.  As  the  user 
moves  the  hot  points,  the  consequent  posture  is 
automatically  predicted/optimized.  That  is,  a  version  of 
the  posture  prediction  algorithm  is  run  with  every  frame, 
using  the  hot  points  as  target  points.  This  version  of 
posture  prediction  currently  has  been  optimized  for 
speed,  and  non-critical  features  have  been  omitted.  The 
algorithm  is  able  to  run  approximately  20  times  per 
second.  This  alleviates  the  need  for  a  user  to  predict 
realistic  postures  as  he/she  positions  the  avatar;  realistic 
postures  are  actually  predicted  automatically. 


Fig.  4.  Inverse  Kinematics 


Fig.  5.  Hot  Points  for  Inverse  Kinematics 


PERFORMANCE  MEASURE  FOR  VISION 

Often,  posture  is  governed  by  a  natural  desire  to  see  the 
target  one  is  striving  to  contact.  This  is  demonstrated 
with  Figs.  6a  and  6b.  Fig.  6a  simply  illustrates  the 
posture  that  is  achieved  when  musculoskeletal 
discomfort  is  minimized.  Although  the  consequent 
posture  is  relatively  comfortable  in  terms  of  joint  stress, 
the  avatar  is  apparently  making  no  effort  to  see  the  target 
point,  hence  the  motivation  for  modeling  the  effects  of 
vision  on  posture. 

This  work  with  vision  is  an  extension  of  the  vision 
function  developed  by  Kim  et  al  (2004).  It  is  based  on 


visual  acuity  and  essentially  has  three  components. 
First,  acuity  is  modeled  as  follows: 

f  _p-7(i -Cose)2 

'acuity  ~  G 

where  6  is  the  angle  between  the  line  of  site  and  the  line 
drawn  directly  from  the  avatar's  head  to  the  target  point. 
In  addition,  a  penalty  function  is  added  such  that  joint 
angles  in  the  neck  avoid  their  extreme  limits.  This 
complex  penalty  function  is  derived  from  the  work  by 
Marler  et  a/ (2005b).  Finally,  tilting  of  the  head  about  the 
axis  that  projects  outward  orthogonal  to  the  avatar's  face, 
is  minimized  such  that  the  corresponding  DOF  gravitate 
towards  a  neutral  position. 

SELF-COLLISION  AVOIDANCE 

A  major  concern  with  human  modeling  is  the  interaction 
between  a  human  and  his/her  environment.  Therefore, 
modeling  self-collision  avoidance  and  collision  avoidance 
with  environmental  obstacles  is  crucial  for  realistic 
posture  prediction.  Obstacle  avoidance  with  multi-body 
systems,  such  as  robots,  has  been  studied  extensively 
(Ozaki  and  Mohri  1986;  Kim  and  Khosla  1992;  Kuffner  et 
al  2002;  and  others).  However,  relatively  little  work  has 
been  done  with  human  self-collision  avoidance  in  the 
context  of  an  optimization-based  approach.  We  have 
implemented  an  approach  with  Santos™,  where  the  self¬ 
collision  avoidance  capability  is  based  on  the  Cartesian 
distances  between  different  body  parts.  This  distance  is 
formulated  as  an  additional  distance  constraint  in  (1).  It 
requires  that  various  end-effectors  on  the  body  do  not 
collide  with  one  another. 


Fig.  6a.  No  Vision 


Fig.  6b.  With  Vision 


Using  spheres  with  various  sizes  to  represent  different 
body  parts,  the  self-collision  avoidance  constraints  are 
expressed  as  follows: 

[xi  (q) — xy  (q)]2  ~{^i+  Rj  )2  (3) 

where  x,.(q)  and  xy.(q)  are  the  Cartesian  coordinates  of 

the  centers  of  spheres  /'  and  j  representing  two  body 
parts,  and  R:  and  4  are  the  radius  of  the  two  spheres, 

respectively.  Distance-squared  is  implemented  here 
instead  of  distance  for  the  convenience  of  numerical 
calculations.  An  example  of  the  effects  of  implementing 
the  constraints  in  (3)  is  shown  in  Fig.  7,  where  a  posture 
is  compared  with  and  without  self-collision  avoidance. 

DYNAMIC  MOTION  PREDICTION 

Dynamic  motion  prediction  is  a  critical  component  for 
virtual  humans.  It  can  generate  joint  torques  and  predict 
joint  injuries.  Hase  and  Yamazaki  (1997)  simulated 
walking  and  rowing  motion  using  muscle  energy 
consumption  formulas  and  human  model,  and  the 
Newton-Euler  method  was  used  to  calculate  the  joint 
torques.  Chaffin  (1997)  developed  two-  and  three- 
dimensional  computerized  human  models  for  simulation 
of  static  strength  during  work.  This  paper  reports  on  the 
development  of  an  optimization-based  dynamic  motion 
prediction  method. 


(  x  =  M(q)q  +  V(q,q)  +  £  J>,g  +  X  VF*  +  K  (q  -  q" ) ) 
where  zt  is  the  generalized  torque  of  joint-/,  Jf  is  the 
transpose  of  the  J  acobian  matrix.  v,(q,q)  is  the  Coriolis 
and  Centrifugal  vector: 


v,M)=ZI  Z  Tr 

k=  1  7=1  _7=max(7,Ar,7) 
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i,k,l  =  \,2,...,n ,  Mik  is  the  (i,k)  element  of  the  mass- 
inertia  matrix  M(q)  such  that 


M„(q)=  Z  Tr 

j=max{i,k) 
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1 
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dqj 
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I,  is  the  inertia  matrix  and  is  written  as  follows: 
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mt  is  the  mass  of  link-/';  (xi,yi,zi )  is  the  location  of  the 

center  of  gravity  of  link-/',  expressed  in  terms  of  t 
coordinate  frame  (local  frame);  and  4,...  are  the 

moments/products-of-inertia  of  link-/' with  respect  to  the  /" 
coordinate  system.  K  is  a  diagonal  matrix  of  elastic 
constants  given  as  follows: 


K  = 


0 


0  kn 

Note  that  all  of  the  variables  are  functions  of  time,  i.e., 
t  =  t (t) ,  Ft  =Ft(f),  and  q  =  q(r) ,  with  k  =  l,2,...,n  .  The 

total  human  metabolic  energy  expenditure  rate  is  the 
sum  of  mechanical  power,  heat  rate,  and  BMR. 

« Zlc-mcol+Z^  |t(o| +b 

7=1  7=1 

B  =  0.685/^  +  29.8  (Joule/second  or  Watt) 
where  W  (Kg)  is  the  body  mass. 


Fig.  7.  Posture  prediction  with  self-collision  avoidance 

The  general  optimization-based  method  is  defined  as 
follows: 

Find:  Control  Points  (PUj  i  =  j  =  \, ...,«) 
to  minimize:  Metabolic  Energy  (EMetabolic) 
subject  to: 

J  oint  limits  (  q\  £  qt  £  qut  i  =  1  ) 

Torque  limits  (?f  £  t,  £  zf  /'  =  1 ) 

Path  constraints  (|| x(q(r))-p(r)  ||^f) 

Equations  of  motion 


The  following  example  is  to  predict  the  motion  of  moving 
a  9  kg  object  on  one  hand  from  a  given  initial  position  to 
a  given  final  position  in  2  seconds  using  the  21-DOF 
Santos™  upper-body  (from  waist  to  right  hand)  model. 
(Fig.  8). 


(a)  At  time  f  =  0  (initial)  (b)  At  time  t  =  0.6 


(c)  Attime  f  =  1.4  (d)  attime  t  =  2  (final) 
Fig.  8.  Santos™  moving  an  object 


The  optimum  joint  kinematic  profiles  are  obtained  at  each 
0.04  of  a  second.  The  joint  torque  values  are  calculated 
and  are  recorded  in  Fig.  9.  For  intuitive  validation  of  the 
above  motion  prediction,  five  human  subjects  were  asked 
to  perform  the  same  task  (with  different  object  weights 
according  to  each  subject’s  maximum  capacity)  and 
demonstrated  similar  motions.  Fig.  10  shows  the  initial 
and  final  posture  of  a  human  subject. 


In  the  initial  posture,  the  spine  joints  for  lateral  bending 
exert  large  torques  for  right  lateral  bending,  whereas  in 
final  posture,  the  spine  joints  for  flexion/extension  exert 
large  torques  for  flexion.  The  large  torque  values  for 
different  motion  segments  will  have  different  effects  on 
the  body  and  may  cause  different  kinds  of  injury  and 
injury  risk  levels  (especially  on  the  low-back).  It  can  be 
observed  that  the  initial  and  final  postures  have  different 
distributions  of  the  torques  at  the  shoulder  and  arm  joints. 
Although  we  have  presented  the  initial  validation  work, 
systematic  validation  using  more  human  subjects  and  the 
statistic  method  is  ongoing. 

MUSCULOSKELETAL  MODELING 

Through  the  understanding  of  simple  mechanical 
principles  and  observing  the  anatomy  of  the  human  body, 
a  model  for  predicting  muscle  forces  to  create  desired 
joint  torques  is  being  developed  for  Santos™.  Developing 
this  model  requires  an  understanding  of  the 
musculoskeletal  system,  the  mechanical  principles  of 
how  torque  is  generated  through  muscle  tension,  and 
optimization  to  solve  the  resulting  indeterminate  problem. 
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Fig.  9.  J  oint torque  profiles 


(a)  Initial  posture  (b)  Final  posture 
Fig.  10.  Visual  demonstration 


The  musculoskeletal  system  is  composed  of  the  skeletal 
system  and  the  skeletal  muscles.  This  system  of  muscles 
and  bones  is  arranged  throughout  the  body  as  numerous 
levers  where  the  bones  act  as  the  levers  and  the  joints 
between  them  serve  as  fulcrums  in  Fig.  11,  which  shows 
some  of  the  major  muscles  of  the  upper  arm  and 
illustrates  how  they  are  attached  to  the  skeletal  system. 


Fig.  11.  Major  muscles  of  the  upper  arm 


Motion  is  then  created  by  the  muscles  pulling  on  the 
bone  to  create  torque  about  a  joint.  Since  muscles  can 
only  act  in  tension,  each  joint  requires  one  set  of  muscles 
to  create  a  positive  torque  and  a  separate  set  to  create 
negative  torque.  In  Fig.  11,  the  elbow  extension  and 
flexion  muscles  are  identified.  Based  on  the  joint  torque 


torque  results  in  the  previous  section,  we  have 
developed  an  optimization-based  method  to  determine 
muscle  forces. 

MECHANICAL  SYSTEM 

Consider  the  simple  one  muscle  model  shown  in  Fig.  12, 
which  represents  the  medial  head  of  the  bicep.  This 
muscle  originates  on  the  scapula  and  inserts  below  the 
elbow  on  the  radius.  For  this  simple  example,  only  the 
torque  created  at  the  elbow  will  be  considered.  When  the 
muscle  is  activated,  it  creates  a  force  on  the  radius.  For 
this  simple  example  it  will  be  assumed  that  this  force 
points  at  the  muscle  origin.  This  force  creates  a  torque 
about  the  elbow  by  acting  on  a  moment  arm  that  lies 
along  the  radius  (Fig.  13). 

This  force  can  be  resolved  into  two  components;  the 
rotational  component  which  is  normal  to  the  moment  arm, 
and  the  stabilizing  component  which  is  parallel  to  the 
moment  arm  (Fig.  14).  The  rotational  component  is  the 
portion  of  the  force  which  actually  contributes  to  joint 
torque.  It  is  important  to  note  that  as  the  joint  angle 
changes,  so  does  the  rotational  component  (Fig.  15).  So, 
the  resulting  torque  is  not  only  a  function  of  the  muscle 
force,  but  also  the  current  joint  angle. For  this  simple  one- 
muscle  example,  determining  the  muscle  force  needed  to 
achieve  a  desired  torque  at  a  specific  joint  angle  is 
simple.  Torque  is  found  with  T  =  rxF ,  where  T  is 
torque,  r  is  a  vector  pointing  from  the  joint  to  the  point 
where  the  muscle  acts  on  the  bone  and  F  is  the  force. 
Unfortunately,  once  more  muscles  are  added  to  this 
problem,  it  is  no  longer  trivial  to  solve.  To  illustrate  this, 
Fig.  16  shows  two  more  muscles  added  to  the  previous 
example. 


Fig.  12.  Medial  Bicep  Brachii 


Fig.  13.  Force  from  Bicep 


Stabilizing 

Component 


Fig.  14.  Muscle  force  components 


Fig.  15.Effect  of  joint  angle 


F  ig.  16  Three  muscles  acting  on  the  lower  arm 


Using  the  same  equation,  torque  is  found  with 
T  =  rbxFb+rbcxFbc+rbrxFbr ,  where  b  denotes  Bicep,  be 
denotes  the  Brachialis,  and  br  denotes  the 
Brachioradialis.  Finding  the  muscle  forces  to  achieve  a 
desired  torque  at  a  specified  angle  becomes  an 
indeterminate  problem  because  there  is  one  equation 
with  three  unknowns  (  and  iv  )■  To  solve  this 

problem,  optimization  must  be  used. 

OPTIMIZATION 

Solving  this  problem  with  optimization  requires  a  cost 
function  which  will  be  minimized  while  obeying  various 
constraints.  There  have  been  many  cost  functions 
proposed,  and  choosing  the  most  appropriate  function  is 
likely  dependent  on  the  situation.  For  example,  quick 
movement  to  avoid  a  collision  will  probably  require  a 
different  cost  function  than  a  movement  that  needs  to  be 
precise  and  smooth.  Many  common  cost  functions  use 
muscle  activation  as  variables.  Muscle  activation  is  a 
value  that  ranges  between  0  and  1  where  0  indicates  no 
muscle  activation  and  1  indicates  a  fully  activated 
muscle.  Muscle  force  can  then  be  defined  as 
F  =  a* FmaY  where  a  is  the  activation  and  Fmav  is  the 
maximum  theoretical  force  that  muscle  can  generate, 
which  is  calculated  by  the  muscle's  physical  cross 
sectional  area.  Torque  can  now  be  represented  as 
T  =  rx(a*Fmax).  To  simplify  the  problem,  only  the  torque 
being  generated  about  the  rotational  axis  of  the  elbow  will 
be  considered.  This  allows  the  torque  equation  to 
become  T  =  r*a*Fm„r  where  r  is  a  constant  and  F  is 
the  rotational  component  of  Fmax.  This  problem  can  now 
be  solved  with  the  following  optimization  formulation 
where  n  is  the  number  of  muscles. 

Find:  al,a2,...,an 

n 

to  minimize:  J  -  y'af 

;= l 

n 

subject  to:  T  =  ^^(a,F;.max)  (4) 

i= 1 

0  <  a,  <  1 

Note  that  by  simply  adding  more  torque  equation 
constraints,  this  formulation  can  be  extended  to  cover 
numerous  joints. 

For  this  example,  ten  muscles  that  are  used  to 
flex/extend  the  forearm  and  pronate/supinate  the  hand 
will  be  considered.  These  muscles  are  the  brachialis, 
brachioradialis,  bicep  long,  bicep  short,  tricep  lateral, 
tricep  long,  tricep  medial,  pronator  teres,  supinator,  and 
anconeus.  Torque  about  the  elbow  and  the  pronating 
axis  of  the  wrist  are  given  and  the  orientation  of  these 
joints  is  known.  The  objective  is  to  find  the  muscle 
activation  level  of  the  ten  muscles  listed  above.  This  will 
be  accomplished  through  optimization  as  described  in 


the  previous  section.  For  this  example,  the  problem  is 
formulated  as  follows: 


Find:  ar..a10 

10 

to  minimize:  J  =  ^]a,2 


subject  to:  TBR  +  TBL  +  TBS  +  Tpt  +  TSUP  =  1 

^br  +  Tbch  +  Tbl  +  ^s.s  +  TTlat  +  ^rio,v  +  ^tm  +  ^pt  Tsw>  +  ^anc  =  1 0 


0  <  a,  <  1 

where  Tx  is  the  torque  generated  by  that  muscle.  Note 
that  for  this  problem,  it  is  given  that  the  torque  about  the 
wrist  is  1  Nm  and  the  torque  about  the  elbow  is  10  Nm.  It 
should  also  be  noted  that  while  all  the  muscle  considered 
contribute  some  torque  to  the  elbow,  only  five  contribute 
torque  to  the  wrist. 


Next  the  torque  generated  by  each  muscle  must  be 
calculated.  These  are  calculated  as 

Ti  =  rtxtatFlm 

where  ris  the  moment  arm,  xis  the  rotational  component 
of  the  unit  vector  representing  the  force  line,  a  is  the 
activation  level,  and  Fmax  is  the  theoretical  maximum 
isotonic  force  which  that  muscle  can  create.  The 
maximum  isotonic  force  can  be  calculated  by  multiplying 
the  maximum  physiological  achievable  muscle  stress 
(600  kPa)  to  the  physiological  cross  sectional  area 
(PSCA)  of  the  muscle  (Berme  et  al,  1987).  The  torque 
equation  can  then  be  represented  as 
=  rixiajPSCAj  *  600 


Note  that  rand  PSCA  are  constants  while  xis  dependent 
on  the  joint  angles  and  a  is  the  design  variable.  Since  the 
joint  angles  are  known,  x  can  be  calculated  from 
geometry  while  rand  PSCA  can  be  found  in  the  literature 
(Maurel  and  Thalman,  1998).  For  example,  TBR  about  the 
elbow  can  be  calculated  as 

Tbr  =  0.0662  *  0.424  *  a  *  4.70  *  600  =  14. 149a 
So,  for  the  given  orientation,  the  constraint  equations  can 
be  calculated  as 

0.002a2  +  1.281a3  +  1.1 12a4  +  (— 0.724)a8  +  0.025a9  =  1 
14.149a,  +8.670a2  +  4.833a3  +  4.228a4  +( — 2.640)a3  +( — 2.617)a8 
+(— 2.636)a7  + 1 ,230a8  +  4.012a9  +  0.01 12a10  =10 


Solving  the  optimization  problem  predicts  the  following 
activation  levels  for  each  muscle:  Bracialis  0.303, 
Brachioradialis  0.186,  Bicep  Long  0.445,  Bicep  Short 
0.384,  Tricep  Lateral  0,  Tricep  Long  0,  Tricep  Medial  0, 
PronatorTeres  0,  Supinator  0.093,  Anconeus  0. 

MUSCLE  STRESS 

Muscle  forces  obtained  in  the  previous  section  are  not 
enough  to  show  all  aspects  of  muscle  injury.  Muscle 
stress  is  one  important  index  for  muscle  injury.  Several 
major  methods  have  been  proposed  for  physically  based 
deformable  modeling  (Gibson  and  Mirtich  1997):  mass¬ 
spring  models,  finite  element  methods,  finite  volume 


models  (Teran  et  al,  2003),  and  other  low-degree 
approximated  continuum  models.  There  are  essentially 
two  issues  in  deformable  body  modeling:  how  to 
parameterize  a  body  and  how  to  bring  in  physical 
behavior  to  the  system. 

Biomechanical  analysis  of  skeletal  muscles  is  an 
important  task  in  the  development  of  a  digital  human 
system.  The  stress  level  within  the  muscle  of  action  is 
particularly  of  interest  to  us.  However,  the  stress  analysis 
with  the  finite  element  method  is  usually  intensive  with 
dense  mesh  as  a  result  of  irregular  shapes  of  muscles. 
Considering  the  variety  of  muscle  shapes  and 
functionality  differences,  simulations  with  over  simplified 
or  idealized  muscle  geometry  and  biomechanics  will  not 
give  satisfactory  results,  especially  when  trying  to  exam  a 
group  of  cooperating  muscles  in  action. 


In  our  work,  a  new  modeling  method  (Zhou  and  Lu  2005) 
is  being  developed  to  overcome  the  difficulties  of 
computation  intensity  of  the  standard  finite  element 
method  in  an  interactive  virtual  environment.  The  method 
is  based  on  the  combination  of  Non-Uniform  Rational  B- 
spline  (NURBS)  geometric  representation  and  the 
Galerkin  finite  element  methods,  which  intuitively  bridges 
the  two  of  them  without  the  meshing  procedure  and 
hence  with  low  degrees  of  freedom. 

We  use  the  NURBS  geometry  to  represent  the  shape  of 
muscles  in  order  to  capture  the  complexity  of  muscle 
geometry.  The  muscle  shape  can  be  extracted  from 
medical  image  data  such  as  the  Visible  Human  Data  Set 
(U.S.  National  Library  of  Medicine).  Our  method  requires 
further  extension  of  NURBS  surface  to  NURBS  solid, 
which  keeps  the  geometric  information  of  surface  while 
extending  to  the  inside.  The  NURBS  solid  representation 
can  be  written  as 

P  (u,v,w)  =  ZZSXA,.  k(u,v,w ) 

/= 0  j=0  k=0 

Using  the  isoparametric  mapping  technique,  the  NURBS 
geometric  shape  functions  are  used  as  the  interpolation 
functions  of  finite  elements.  Therefore,  the  displacement 
mapping  can  be  written  as 

n  m  i 

d(u,v,w)  =  vvv,.  R  .  k(u,v,  w) 

i= 0  7=0  k=0 

where  is  the  displacement  of  the  control  point  p  . 

Based  on  the  above-established  mapping,  we  can 
calculate  the  deformation  gradient  and  other  needed 
quantities  for  the  finite  element  computation.  Because 
muscles  are  general  under  very  large  deformation  during 
human  motions,  a  fully  nonlinear  formulation  for  the 
analysis  of  large-strain  motion  of  muscle  is  used. 

CONSTITUTIVE  MODEL  OF  SKELETAL  MUSCLE 

The  internal  structure  of  skeletal  muscle  is  comprised  of 
soft  tissues  with  various  material  properties.  More 
important,  skeletal  muscles  are  active  and  the  properties 


of  the  active  tissue  are  altered  upon  activation.  In  our 
work,  the  muscle  model  is  described  by  a  hyperelastic, 
anisotropic  constitutive  equation  that  includes  an  active 
component. 


and  a  swing  dynamic  prediction  model  is  being 
developed.  Collision  avoidance  is  being  extended  from 
self  collision-free  to  environment  collision  avoidance. 
Finally,  a  user-  friendly  interface  is  being  developed. 


The  hyperelastic  strain  energy  function  is  assumed  to  be 
the  sum  of  two  parts:  energy  function  of  the  passive 
ground  substance  and  active  fibrous  structure.  The 
passive  ground  substance  consists  of  connective  tissue, 
water,  etc,  and  is  modeled  as  isotropic  neo-Hookean 
material.  The  muscle  fiber  strain  energy  is  assumed  to  be 
a  function  of  the  muscle  fiber  stretch  and  the  muscle 
activation  level.  The  muscle  fiber  directions  are  generally 
distributed  with  different  patterns  among  different 
muscles.  In  our  work,  we  used  the  isocurves  of  the 
NURBS  solid  to  model  the  muscle  fibers  and  its 
directions. 

RESULTS  OF  MUSCLE  STRESS 

The  motion  analysis  provides  at  discrete  points  of  time 
the  joint  positions  and  rotations  as  well  as  the  activation 
level  in  individual  muscles.  From  the  position  of  bones, 
the  displacement  boundary  conditions  for  the  muscles 
can  be  known.  And  the  activations  of  muscles  are  fed 
into  the  constitutive  model.  The  stress  analysis  is  then 
conducted  using  the  NURBS-FEM  developed.  As  shown 
in  figure  1,  both  passive  and  active  motions  of  muscles 
are  simulated,  and  the  corresponding  muscle 
deformations  and  stresses  can  be  obtained. 

Using  the  pseudo-static  simulation  of  muscle  deformation 
at  discrete  points  of  time,  the  muscle  force-extension 
curves  are  shown  in  Fig.  17  with  different  length 
changes.  The  curves  (Fig.  18)  show  a  behavior  that 
qualitatively  matches  what  is  reported  in  the  literature 
(Zajac  1989).  It  can  be  seen  that,  with  the  same  length, 
the  active  muscle  generates  more  force  than  the  passive 
one,  and  the  largest  active  muscle  force  generated 
happens  at  the  initial  rest  length  of  muscle. 


SUMMARY  AND  DISCUSSION 

This  paper  presents  the  new  capabilities  of  the  virtual 
human  (Santos™).  These  new  capabilities  include  hand 
modeling,  advances  in  posture  prediction  (new  vision 
objective  function,  multiple  end-effectors,  advanced 
inverse  kinematics  (IK),  self-collision  avoidance), 
dynamic  motion  prediction,  and  musculoskeletal 
modeling  (optimization-based  method  to  determine 
muscle  forces  and  NURBS-FEM  approach  for  muscle 
stress). 

In  addition  to  the  above-mention  newly  developed 
capabilities,  exciting  research  is  still  on-going  and 
concerns  more  features  for  Santos™.  Particle  and  FEM 
methods  are  implemented  to  clothing  modeling. 
Ingress/egress  is  studied  using  motion  capture  system 
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Fig.  17.  Stress  analysis  of  muscle,  (a)  initial  rest 
geometry  (b)  active  contraction  (c)  large  passive 
deformation. 


Fig.  18.  Force-extension  curve 
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